Non-equilibrium GW Approach to Quantum Transport in Nano-scale Contacts 
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Correlation effects within the GW approximation have been incorporated into the Keldysh non- 
equilibrium transport formalism. We show that GW describes the Kondo effect and the zero- 
temperature transport properties of the Anderson model fairly well. Combining the GW scheme with 
density functional theory and a Wannier function basis set, we illustrate the impact of correlations 
by computing the I-V characteristics of a hydrogen molecule between two Pt chains. Our results 
indicate that self-consistency is fundamental for the calculated currents, but that it tends to wash 
out satellite structures in the spectral function. 
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PACS numbers: 72.10.Bg,73.23.-b,73.63.Rt 



Electronic correlations are responsible for important 
transport phenomena such as Coulomb blockade and 
Kondo effects [lj, yet its significance for transport in 
nano-scale structures is not well understood nor has it 
been systematically studied. At present, the most pop- 
ular approach to ab intio simulations of transport in 
nanocontacts combines a non-equilibrium Green's func- 
tion formalism with the single-particle Kohn-Sham (KS) 
scheme of density functional theory (DFT). This ap- 
proach works well for some systems @, Q , but in other 
cases it fails to reproduce experimental data Q indicat- 
ing the need for computational transport schemes beyond 
the DFT level 

A reliable description of transport through a molecu- 
lar junction requires first of all a reliable description of 
the electronic structure of the molecule itself, i.e. its elec- 
tron addition- and removal energies. It is well known that 
the GW self-energy method yields q uasip article proper- 
ties of molecules d, Q and solids [l(J in good agree- 
ment with experiment improving drastically the DFT 
band structures. In view of this it seems tempting to 
extend the use of the GW approximation to transport 
calculations. It is clear, however, that this should not 
be implemented by shifting the molecular energy levels 
to their GW positions prior to coupling. The reason is, 
that when a confined interacting system is connected to 
external (non-interacting) leads, the electrons in the con- 
fined region become correlated with those in the leads. 
To capture these correlations, which are the origin of im- 
portant many-body phenomena such as the Kondo effect, 
it is crucial that the self-energy be evaluated in the pres- 
ence of coupling to the leads. 

Traditionally, correlation effects in transport have been 
studied on the basis of the Anderson- and Kondo mod- 
els by a variety of numerical and analytical techniques. 
Many of these techniques are, however, quite specific to 
the considered models and lack the generality needed to 
be combined with first-principles methods. 

In this paper we combine the GW approximation with 
the non-equilibrium Keldysh formalism to obtain a prac- 



tical scheme for correlated quantum transport. We study 
the Anderson model out of equilibrium, and calculate the 
I-V characteristics of a molecular hydrogen contact using 
a Wannier function (WF) basis set. In both applications 
we emphasize the difference between self-consistent and 
non self-consistent evaluations of the GW self-energy. 

As a general model of a quantum conductor we con- 
sider a central region (C) connected to left (L) and right 
(R) leads. The leads are kept at chemical potentials 
fiL and fiR, respectively. We construct the matrix 
hj = (cf>i\h s \(t>j), where h s = + v h + v xc + v ext 

is the KS Hamiltonian of the combined L-C-R system 
in equilibrium, and {<f>i} is a corresponding set of max- 
imally localized, partly occupied WFs [12 |. Assuming 
that correlation effects as well as charge redistributions 
induced by the bias voltage are significant only inside 
C, we describe the leads and the couplingto the central 
region by h. Using standard methods [1 31 ] we evaluate 
the coupling self-energies, E a (w) = hc a 9ua{u)h a c, 
where g aa is the GF of the uncoupled lead a = L,R. 
The interactions inside the central region are de- 
scribed by Vi nt = \Y, l jkl£C,aa' V V,kl4a C )<y' Cl °' Ck °i 

with the Coulomb matrix elements Vij u = 
J/drdr'^(r)*^(r');^(r)^(r')/|r - r'|. In prac- 
tice we use an effective interaction which only involves 
a subset of the Coulomb matrix elements, see later. 



We include matrix elements of the form V% 



and 

in the calculation of the correlation part of the 
GW self-energy while also terms of the form Vujj 
and Vu y ij are included in the Hartree and exchange 
self-energies [l4| • By allowing the effective interaction to 
be spin-dependent we avoid the self-interaction normally 
present in the GW correlation. 

We calculate the retarded and lesser Green's functions 
of the central region from [HI, [l6| 



G" 

G< 



[w + irj - hoc + (v X c)cc - &v h - E[ ot ] 1 (1) 
G r Z< t G a + 2r ] ifcG r G a , (2) 



where 77 is a positive infinitesimal and Utot = El + S/j + E 
is the sum of the coupling self-energies and the exchange- 
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correlation part of the interaction self-energy. The term 
Av h = E£[G] — E£[Gp FT ], is the correction to the equi- 
librium (DFT) Hartree potential introduced by the cor- 
relations and the finite bias. The DFT xc-potential is 
subtracted from hcc to avoid double counting. In the 
last term of Eq. ([2|), fc denotes the initial Fermi-Dirac 
distribution of the central region before coupling to the 
leads. Notice that this term can become significant at 
energies where ImE[ ot is comparable to r\ While 
this is not expected to occur in the interacting case (due 
to life-time broadening by E), it will happen in the non- 
interacting case whenever bound states are present in G. 
We note, that in Eqs. (P) and ([2]) we have specialized 
to the long time limit where we assume that the system 
reaches a steady state in which the Green's functions de- 
pend only on the time difference t = t% — 1\ and thus can 
be represented by a single time/frequency variable. (l7j 

The symmetrized current, / = (Il-\-I r )/2, whe re I ta r) 
is the current in the left (right) lead, is given by [l8j 

where ^l(r) — *Pl(_r) ~ ^i(ii)] ^ s t ne coupling strength 
and the trace is taken over basis functions in the central 
region. 

Within the GW approximation E is written as a prod- 
uct of the Green's function, G, and the screened in- 
teraction, W, calculated in the random phase approx- 
imation (RPA). Out of equilibrium this holds true on 
the Keldysh contour, and the relevant equations in real 
time follow from the Langreth conversion rules [15j | . Ab- 
sorbing the spin index into the orbital index we define 
the effective interaction V c g = \ ^2 i jVijc\c S jCjCi 1 where 



Vi 



y ,%} 



Kff resembles the real space in- 
teraction with the important difference that Vij is spin- 
dependent and Vu = 0. Self-interaction is thus automat- 
ically excluded to all orders in V e g. The retarded and 
lesser GW self-energies become (on the time axis), 

Er.(t) = iG^(t)W^(t) + iG^(t)W^(t) (3) 
E<(i) = iG£.(t)W£(t). (4) 

The screened interaction is given by (in frequency space) , 

W r {uj) = V[I - P r '{uW]- 1 (5) 
W </> (uj) = W r {Lu)P </> {Lu)W a {Lu), (6) 

where all quantities are matrices in the central region 
indices. Finally the irreducible polarization becomes 



p </> 

ij 



Pli(t) 
(t) 



^(t)G<(-t) - iG<(i)G£(-t) (7) 

(8) 



-iG<'>(t)G>/<{-t). 



In principle the GW method implies a self-consistent 
problem, i.e. the G obtained from Eqs. (fT]) and © should 
equal the G used to produce E. In this way all impor- 
tant conservation laws are fulfilled 191. However, due 



to the large computational demands GW band struc- 
ture calculations usually apply a, quite successful, non 
self-consistent Go Wo approach using on the Go obtained 
from a KS calculation [lC]]- We mention here that the 
non-equilibrium GW approximation has previously been 
used in the study of semiconductors in high-intensity 
laser fields [HHof. 

We represent all quantities on a uniform real 
time / frequency grid using the Fast Fourier transform to 
switch between the two representations to avoid time con- 
suming convolutions. The calculation of the WFs, the 
matrix elements hij , and the coupling self-energies E Q is 
described in detail in Ref. [l3|. 

Before turning to the Wannier-GVF calculations we ap- 
ply the method to the Anderson impurity model. Al- 
though GW is not expected to be accurate for strongly 
correlated systems, the simplicity of the Anderson model 
makes it ideal for illustrating the properties and limi- 
tations of the non-equilibrium GW approximation. We 
thus consider a central site of energy e c coupled to onc- 
dimcnsional leads with on-site energy Eq — and nearest- 
neighbor hopping to = 10. The hopping to the central 
site is t c = 1.8 giving T = T L {0) = T R (0) w 0.65. Double 
occupation of the central site costs a charging energy of 
C/ = 4, and the spin-dependent interaction entering Eq. 
([5|) is V aC r' = U(l — 8 aa i). We assume half-filled bands, 
i.e. Ep — So, and measure all energies relative to Ep- 

We consider three different approximations: (i) Self- 
consistent Hartree, which is equivalent to self-consistent 
Hartree-Fock (HF) since the exchange term vanishes, (ii) 
Go Wo with the self-consistent HF Green's funtion as Go- 
(iii) Fully self-consistent GW. All calculations are non- 
magnetic, i.e. Gf = Gj = G. 




FIG. 1: Spectral function at the central site, e c = —3, of the 
Anderson model calculated in three different ways (see text). 
The curves have been vertically offset for clarity. Inset: real 
and imaginary parts of the GW self-energy together with the 
line ui — ehf- The steep shape of ReS r around Ef pins the 
position of the quasi-particle peak. 

In Fig. [T] we show the equilibrium spectral function at 
the central site, ImG r (w), for e c = —3. The HF solu- 
tion shows a single peak at £hf = £c + U{h a ) with a 
full width at half maximum given by 2r. This behavior 
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is representative for any mean-field description, includ- 
ing the KS scheme. The inclusion of dynamic correla- 
tions leads to qualitative changes in the spectral peak 
which moves close to the chemical potential [2f[ and nar- 
rows down from 2I 1 to 0.63 (0.28) in the case of Go Wo 
(GW). This change is a signature of the Kondo effect: 
For r — U < e c < — T (the so-called Kondo regime), the 
correlated groundstate is a singlet with a finite ampli- 
tude for the central site being empty. At T = this 
leads to the formation of a spectral peak at the chemical 
potential with a width given by the Kondo temperature, 
T K = i(2r[/) 1 / 2 exp[ 7 r£ c (e c + U)/2TU]. For our choice 
of parameters Tk = 0.19, which is in fair agreement with 
the GW result, and about three times smaller than the 
GoW) result. 

From the inset of Fig. [1] it can be seen that the Kondo 
peak gets pinned to E F due to the steep shape of Re£ r 
in this region, and that its reduced width, as compared 
to 2r, is a consequence of the steep drop in ImE r away 
from the point ImS r ((jj = Ef) = 0. The atomic levels 
which should be seen at s c and e c + U appear as shoul- 
ders in the Go Wo spectrum, but for larger values of U/T 
they become more pronounced as satellites (side-bands) 
of the main quasiparticle peak although their positions 
are somewhat shifted towards Ep. In contrast the self- 
consistent GW fails to capture the side-bands. These 
findings agree well with previous results obtained with 
the fluctuation-exchange approximation [2^ |. and with 
GW studies of the homogeneous electron gas (23|. We 
mention that in self-consistent second-order perturbation 
theory, the pinning of the main spectral peak to Ep is 
less pronounced than in GW and its width is significantly 
overestimated, showing as expected that the higher order 
RPA diagrams enhance the strong correlation features. 

In Fig. [2] we show the zero-temperature differential con- 
ductance under a symmetric bias, (J.l/r — iV/2, as a 
function of e c . For V = there is only little difference 
between the HF and GW results which both shows a 
broad conductance peak reaching the unitary limit at the 
symmetric point e c — —U/2. The physical origin of the 
conductance trace is, however, very different in the two 
cases: While the HF result is produced by off-resonant 
transport through a broad spectral peak moving rigidly 
through the Fermi level, the GW result is due to trans- 
port through a narrow Kondo peak which is always on 
resonance (for e c in the Kondo regime). This difference 
is brought out clearly as V is increased: for V <C T the 
bias has little effect on the HF conductance while the 
GW conductance drops dramatically already at biases 
comparable to Tk due to suppression of the Kondo reso- 
nance. We note that we do not observe a splitting of the 
GW Kondo resonance at finite V The side peaks 

in the dl/dV correspond to the central level crossing the 
chemical potentials. The GoW) conductance is markedly 
different from the HF and GW results. At low bias there 
are unphysical dips in the conductance curve, and as V is 
raised the dl /dV becomes even more unreasonable show- 
ing strong negative differential conductance (for this rea- 



son we have omitted the V = 4 curve). This unphysical 
behavior of the GoW) is a result of its non-conserving na- 
ture, and it underlines the necessity of using a conserving 
approximation such as the self-consistent GW for trans- 
port calculations. Within GoW) the violation of current 
conservation, (Jj, —Ir)/ I, increases with V reaching 20% 
at V = 0.8 for certain values of e^, while it is always neg- 
ligible in our GW calculations [25j |. 

As an illustration of the Wannier-GW scheme we 
consider a molecular hydrogen bridge between infinite 
atomic Pt chains, see inset of Fig. [3] Experimentally, 
the conductance of the hydrogen contact is found to 
be close to the conductance quantum, 2e 2 /h, and this 
value has been reproduced by DFT calculations. [|| Us- 
ing the plane- wave pseudopotential code Dacapo, we per- 
form DFT calculations for an infinite Pt wire as well as 
a supercell containing the hydrogen molecule with six Pt 
atoms on each side. The WFs and KS Hamiltonian of the 
full system in equilibrium are obtained by combining the 
two calculations, see Ref. [HI for more details. For the 
transport calculation we use a central region consisting of 
the two s-like WFs of the hydrogen. Due to the smallness 
of the matrix elements coupling the two Pt chains across 
the molecule this suffices to converge the elastic trans- 
mission function. Ultimately, the dependence of the GW 
results on the size of G as well as on the basis set should 
also be checked. (26j . 

In the upper panel of Fig. [3] we show the local density of 
states (LDOS) at one of the two H orbitals as calculated 
within DFT using the PW91 xc- functional, as well as 
self-consistent HF (in the central region) . In DFT the H2 
bonding state is a bound state at —7.0 eV relative to Ef, 
while the anti-bonding state lies at 0.4 eV and is strongly 
broadened by coupling to the Pt. Moving from DFT to 
HF the bonding state is shifted down by ~ 8 eV because 
for occupied states the exchange potential is more nega- 
tive than the DFT xc-potential. The same effect tends 
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FIG. 2: (color online). Differential conductance in the An- 
derson impurity model as a function of the central site en- 
ergy, e c , for different applied biases. The negative differential 
conductance seen in the middle panel is an artifact of the 
non-conserving nature of the Go Wo approximation. 
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E-E F (eV) 

FIG. 3: (color online). Local density of states at one of the 
H orbitals of the Pt-H-H-Pt contact shown in the inset. 
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FIG. 4: (color online). I-V and di jdV for the hydrogen con- 
tact as calculated in DFT(PW91) and self-consistent GW. V 
is the source-drain bias voltage. 

to drive the half-filled anti-bonding state down but in 
this case the resulting increase in the Hartree potential 
(about 4 eV) stops it just below Ep. 

In the lower panel of Fig. [3] we show the LDOS calcu- 
lated in GW as well as GqWq starting from either DFT 
or HF, i.e. Go is either Gdft or Ghf- The large de- 
viation between the two GqWq results is not surprising 
given the large difference between Gdft and Ghf • Focus- 
ing on the bonding state, the GqWq quasiparticle (QP) 
energies lie at —26 eV and — 11 eV for DFT and HF, 
respectively. Two effects are responsible for this differ- 



ence: We have eqp(Go) = £hf(Go) + s corr (Go), where 
e corr is determined by the intersection of the line s — £hf 
with the real part of the correlation self-energy, ReT, r corr . 
Now, £hf(Ghf) is already ~ 5 eV larger than £hf(Gdft) 
due, mainly, to the mentioned increase in Hartree en- 
ergy. Secondly, it turns out that T, r corr (GHF) is roughly 
S^ orr (GDFT) shifted down by ~ 9 eV (note that this 
corresponds to the difference between the delta peaks in 
Ghf and Gdft), leading to e CO rr(GDFT) ~ —6 eV and 
£corr (Ghf) ~ 5 eV [HI- We are aware that part of this 
difference could be due to the limited size of the basis. 

We notice, that the LDOS results of Fig. [3] can be 
largely reproduced by including only the second-order 
GW diagram in the self-energy. Thus the higher-order 
RPA diagrams are less important in this case. 

The linear-response conductance has been calculated 
by applying a small bias of 10 mV. All the self-consistent 
calculations, i.e. DFT, HF, and GW, yield a conduc- 
tance within 10% of the experimental value of 2e 2 /h. 
The same holds for GqWo(HF), however, this is some- 
what arbitrary as the Go Wo (DFT) conductance is only 
0.4(2e 2 //i). In Fig. 2] we show the fully self-consistent 
I-V characteristics for DFT and GW. The DFT con- 
ductance is nearly constant over the bias range (like the 
HF, not shown). In contrast the GW conductance falls 
off at higher bias due to incoherent scattering described 
by JmS corr . Since ImY. corr (Ep) vanishes in equilibrium, 
the finite-bias conductance suppression is a direct result 
of the non-equilibrium treatment of correlations. 

In conclusion, we have presented a non-equilibrium 
GW implementation that can be combined with DFT 
and a localized basis set to model correlated electron 
transport in nanostructures. Results for the Anderson 
model and a Pt-H-H-Pt nano-contact indicate that self- 
consistency is crucial for GW transport calculations. 
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